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Abstract 

Recent observational surveys of trans-neptunian binary (TNB) systems have dramatically increased the number of 
known mutual orbits. Our Kozai Cycle Tidal Friction (KCTF) simulations of synthetic trans-neptunian binaries show 
that tidal dissipation in these systems can completely reshape their original orbits. Specifically, solar torques should 
have dramatically accelerated the semimajor axis decay and circularization timescales of primordial (or recently ex- 
cited) TNBs. As a result, our initially random distribution of TNBs in our simulations evolved to have a large popu- 
lation of tight circular orbits. This tight circular population appears for a range of TNO physical properties, though a 
strong gravitational quadrupole can prevent some from fully circularizing. We introduce a stability parameter to pre- 
dict the effectiveness of KCTF on a TNB orbit, and show that a number of known TNBs must have a large gravitational 
quadrupole to be stable. 
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1. Motivation 

Trans-neptunian binary systems (TNBs) constitute 
at least 10% of th e objec ts between 30 and 70 AU 
dStephens and Noll 120061) . and up to 30% o f the 



2008b). As 



Cold Classical Kuiper Belt jNoll et al 
of spring 2012, 72 TNBs have been reported in 
the literature, with full mutual orbits having been 
reported for 18 objects, partial orbits with am- 



biguous orbi ts for 3 mor e (e.g. iNoll et al.L 12 008a: 



Grun dy et al.L 120091 1201 U iParker et al.L 1201 ll and 



the list at http://www2.lowell.edu/users/grundy/tnbs). 
These observations show that the majority of detected 
TNB systems have a separation of less than 2% of the 
Hill Radius {rnm), defined as: 



rmii = aheiioi ^ - eheiio) 



,4 



3M, 



(1) 



S un 



where a/, e / 1£ , and eheiio are the semimajor axis and ec- 
centricity of the heliocentric orbit. Even more striking 
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is the very small fraction of TNB systems which are 
widely-separated (>10% a/rnm), despite their being 
easier to detect. This implies that TNBs are generally 
in very close mutual orbits, and the fraction of orbits 
that are very close has only increased with better 
detection methods. In addition, most known TNBs 



are of almost equal brightness dNoll et al 1 l2008ah 
implying near-equal masses. 



Several formation methods have been proposed to 
create TNBs, though none as yet can fully describe 
the observed population, nor account for any post- 
formation orbital evolution. Large impacts are an 
obvious contender for formation, but tend to produce 
smaller satellites (and thu s less e qual mass ratios) 
than are observed. Indeed, Canurj (2005) showed the 
Charon-forming impact required a very slow relative 



velocity (v,-, 



0.7 km/s), and even then only 



allowed a mass ratio of approximately 1 0:1. Dynamical 



captures can also prod uce TNBs (e.g. Goldreic h et al. 



2002; ILee et all 120071) . These methods do favor near- 
equal mass ratios (as it provides a deeper gravity well 
per size of the primary object), but have great preference 
for pro ducing wide ( >5% rum) binaries on eccentric 
orbits. Funato et al. I d2004 combines a small impact 



Preprint submitted to Icarus 



March 8, 2013 



and dynamical capture to efficiently produce T NBs, 
but on ly at very high eccentricities. iNesvornv et al. 



(2010) shows that binaries formed by gravitational 
collapse also tend to have near-equal mass ratios, but 
again have wide, moderately eccent ric orbits. The 
unbinding of binaries by impacts (Petit and Mousis, 
2004) or Neptune encounters (IParker and Kavelaars 
2010) would reduce the number of wide TNBs, but 
would not correspondingly increase the number of tight 
systems. The deficit of these wide systems, and the 
abundance of tight ones, therefore hints at the existence 
of some non-disruptive post-formation processing of 
TNB mutual orbits. 

In this paper, we pr opose that Kozai Cycle Tidal Fric- 
tion (KCTF, after lEggleton and Kisseleva-Eggletonl 
2006) may be the method by which these orbits were 
tightened and circularized. We will show through 
several sets of Monte Carlo simulations that KCTF can 
transform a large fraction of primordial TNB systems 
into very close and circular orbits. In addition, we show 
that the tidal evolution this implies means that Kozai 
cycles are very inefficient at destroying TNB systems. 
We also show how KCTF is influenced by the physical 
properties of the TNB system, such as tidal Q and ki, 
density, J2, rotation rate, and mass ratio. 

Figure 1 here 



2. KCTF Model 

In order to understand how TNB orbits may have 
evolved since they were formed, we created a numerical 
Kozai Cycle and Tidal Friction model. Kozai Cycles in 
this context are periodic oscillations in eccentricity and 
inclination of the TNB mutual orbit caused by solar 
torques. For this paper, the outer orbit is the heliocentric 
orbit of the binary's barycenter, and the inner orbit is 
the mutual orbit of the binary pair. These oscillations 
preserve the orbit's semimajor axis and the quantity 

? 2 , where / is the inclination of the 



cos I x Jl 

V in 

mutual orbit with respect to the heliocentric orbit and 
e ln is the eccentricity of the i nner o rbit. This process 
was first described in iKozail (119621) . in the context of 
perturbations by Jupiter on asteroid orbits. Without any 
tidal or quadrupole effects, these oscillations would 
vary eccentricity periodically with a period between 
approximately 2 ka and 2 Ma. IKozail (1 19621) showed 
that for an initially circular orbit, the minimum inclina- 
tion for oscillation to initiate is ±39.2°. However, this 
limiting inclination becomes much lower at non-zero 



initial eccentricities. Thus, the only mutual orbits that 
could be excluded from this effect are those which form 
with both low initial eccentricity and low inclination 
relative to their heliocentric orbit. Kozai cycles have 
been suggest ed as a method of evolv ing the mutual or- 
bits of TNBs ( Perets and Naoz , 20091) . as well as binary 
and triple Near-Earth asteroids, whic h may have Koza i 
cycles short enough to be observable (IFang et al.L 1201 lb . 



A significant consequence of these Kozai oscilla- 
tions is that the eccentricity of the mutual orbit can 
become very high, especially if the initial orbit has a low 
eccentricity but high inclination (or vice versa). Since 
the tidal dissipation rate for these objects is chiefly a 
function of their mu tual separation at periapse (see 
Equa tions 5 to 7 in Eggleton and Kiseleva-Eggletonl 
2001), a minor increase in eccentricity can have a 
major effect on the amount of tidal dissipation. This 
is important, as tidal models that assume near-zero 
eccentricity would produce much slower tidal evolution 
than is realistic for an orbit with high eccentricity due to 
solar Kozai effects. Mutual orbits with Kozai-pumped 
eccentricities can therefore decay due to tidal friction 
much faster than their initial state would imply; see 
Figure Q] for an example. The strength of Kozai-driven 
tidal decay is inversely proportional to the magnitude 
of the binary orbit's angular momentum as projected on 
the axis of the heliocentric orbit's angular momentum 
vector. Also, because this projected angular momentum 
is perpendicular to the solar-driven precession of the 
system, it can be completely determined from the in- 
stantaneous orbit without knowledge of the precession. 
Here, we normalize this quantity to a percent of the Hill 
radius of the system: 



H' = cosI A \a in (l - 4,)— — 
rmu 



(2) 



Where a m is the semimajor axis of the inner orbit. We 
find this to be a particularly useful normalization, as 
values of H' smaller than one will experience strong 
body tides over the course of a Kozai cycle, while 
larger values generally will not (depending on their 
physical properties). Also, since we do not evolve the 
binary's heliocentric orbit in our simulations, lOQ/rnni 
is a constant normalization parameter. H' is effectively 
Tisserand's Parameter for a three-body system with 
only quadrupole perturbations, which is appropriate 
here because all known TNBs have a/, e / !0 »a,„ by at 
least five orders of magnitude. 

Since H' is a much stronger function of the orbit's ori- 
entation than separation, even very wide binaries can be 



2 



affected by KCTF if their inclination (or eccentricity) is 
high enough. This KCTF process has been previously 



identified as significant for TNBs by Perets and Naoz 



(2009), but only demonstrated for the Orcus-Vanth 

J 1 1 71. 

dwarf planet system (Ragozzine, 2009). We used a sim- 
ilar m odel based on Eggleton and Kisele va-Eggleton 
J200l[ herein EKE01) to iRagozzind fe009b. which 
is described below. This model directly evolves the 
mutual orbital elements and spin vectors of the binary 
while holding the heliocentric orbit constant. We 
did not include any dynamical effects from objects 
external to the binary other than the Sun. The general 
eq uations of the EKE01 KCTF m odel are summarized 
in iFabrvckv and Trem aine d2007h : below we describe 
the modifications and additions we used. These consist 
of our estimation of frictional timescale, quadmpole 
gravity terms, and integration methods. 

2.1. Frictional Timescale 

Since the EKE01 model was developed for binary 
stars and giant planets, we needed to modify the terms 
relating to the physical characteristics of the objects. 
One benefit of the EKE01 method is that all these terms 
are condensed into a single frictional timescale for each 
object. This timescale, however, is also a function of the 
mutual orbit's semimajor axis and is thus time-varying. 
We therefore reformulated the frictional timescale in a 
way that is computationally more useful. 

In the EKE01 model, the behavior of the objects' 
body tides is determined by the second tidal Love 
number (ki = k<i) and the tidal dissipative function (Q) 
for each object. The Love number is highly dependent 
on both the composition of the object and whether it 
is physically a solid object or a rubble pile. For half 
of our simulations, we assu med the obj ects were solid 
homogeneous elastic bodies (Burns, 1977): 



, _3/ l9fi r R \~ 

klmu - 2 \y + 2GMp j 



(3) 



We took the rigidity of the objects, /j. r , to be 
4 X 10 9 N/m 2 , using the value for icy bodies from 



Gladman et al 



(1996). For the other half of the simu- 
lations, we assumed t he objects were rubble pile s, using 
the approximation of iGoldreich and Sari (200§>: 



R 



<L,rubble - 



10 5 km 



(4) 



In addition, we assumed a constant value for Q, the in- 
verse of the average fr action of tidal energy lost to heat 
per radian of the orbit (Goldreich and Soter 1966). We 



can then find the tidal timescales as a function of the bi- 
nary's orbit (a,„,n,„), ki, and Q. The viscosity (fy) and 
frictional fa) timescales for the primary obje ct as were 
formulated in lFabrvckv and Trema ine (2007, Equations 
A9 and A10) as: 



3 (1 +k u ) 2 Qm in R\ 
2 k LA GM\ 



tv.l 

~~9 



Rij (Mi + M 2 )M 2 



(l+k LA ) 



-2 



(5) 
(6) 



The timescales for the secondary object are the same 
equations, but switch the subscripts 1 and 2. Combin- 
ing these two equations allows a,„ to be separated out, 
thus reducing the number of calculations required per 
iteration: 



i e, m x 



R: 



6 k LJ M 2 VG(Mi + M 2 ) 



x a. 



13/2 



(7) 



Note the leading fac t or of 1 /6 is erroneously listed 
as 2/3 in Rag ozzine d2009l) . leading to longer fric- 
tional timescales, and thus slower orbital decay. In 
addition, this timescale is for an object in a perfectly 
circular orbit with its rotation synchronized to the orbit. 
Equations 5 and 6 in EKE01 combine these factors to 
account for the eccentricity of the orbit and the rotation 
of the objects. Since the closest separation of the two 
objects is the key driver for tidal evolution, the effective 
timescale is very sensitive to eccentricity. 

This frictional timescale gives the approximate 
rates of evolution for a near-circular orbit, but the actual 
rates are strongly dependant on the orbit's eccentricity. 
To illustrate this point, consider the special case of 
a system where the objects have equal values of t F , 
have e » 0, and their rotation is synchronized to the 
orbit, the rates of change for the semimajor axis and 
eccentricity can be approximated as: 



da -a 

~dt ~ t F (l-e) 15 ' 2 

de -1 

dt ~ t F (l-e) 13 ' 2 



(8) 
(9) 



While the eccentricity is e » 0, the semimajor axis de- 
cay can then be estimated as: 



a 



ff,o(l - e ) 



15/2 



f + f F ,o(l-eo) 15/2 



(10) 



Where ao, eo, and tf o are the initial semimajor axis, 
eccentricity, and frictional timescale. As an example, 
consider an equal-mass TNB system with objects 
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having a radius of 100 km, 2=100, rubble-pile ki, 
and density of 1.0 g/cm 3 . At a semimajor axis of 
10 4 km (100 radii), the objects would have values of 
tp » 7 X 10 12 years. If the initial eccentricity were 0.5, 
the orbit would only have decayed to 9998 km after one 
million years. At an eccentricity of 0.8, the semimajor 
axis would decay to 8512 km in that time. And at an 
eccentricity of 0.9, the orbit would shrink to just 306.5 
km (3 radii) after a million years. By the time the 
system reaches an orbit this tight, tf has decreased to 
less than 1000 years, allowing for rapid circularization. 
Clearly then, only a brief excursion to high eccentricity 
is needed to start a feedback loop of tidal decay to a 
tight circular orbit. In KCTF, those excursions happen 
when Kozai cycles pump up the eccentricity. 

The system in Figure Q] shows this process in ac- 
tion. The two objects are equal mass rubble-piles with 
Q=10 and radii of 42 km. The initial orbit is at 9.7% 
of rjiiiu e=0.99, and inclination of 99°. However, the 
system's orientation puts it initially past the peak of 
its Kozai cycle, and so the orbit starts to become more 
circular and less inclined with little tidal evolution. 
After reaching a minimum of <?=0.46 at 3200 years, the 
eccentricity then grows again. By the time it reaches 
e=0.98, the semimajor axis is shrinking at a rate of 2 
km/year and 45 km/year at e=0.983. The shrunken 
semimajor axis reduces tf , so the peak decay rate is 
150 km/year at 6550 years since the start, e=0.985, and 
a-6.1% of rfiui. This corresponds to a periapse distance 
of 212 km, or 5 radii apart. The orbit then gradually 
decays down to a tight circular orbit at 0.15% rum or 
8.1 radii. 

Most KCTF systems do not evolve as rapidly as 
this, and will often have several eccentricity peaks 
with a small amount of decay until the periapse be- 
comes close enough for rapid decay to occur. Still, a 
large range of orientations can produce strong KCTF 
evolution. To measure this effect, we reran a number 
of the simulations described below without any solar 
perturbations. In all these cases, the number of systems 
which circularized was reduced by at least a factor 
of 3.5 compared to the full KCTF model. These 
systems all oriented randomly on the sphere of the 
sky (corresponding to observed systems), so the joint 
process of KCTF can cause fast tidal evolution in a 
large variety of TNB mutual orbits. 

Figure 2 here 



2.2. Quadrupole Gravity 



Ragozzine 009) panded this model by adding the 
capacity for the objects to have a permanent quadrupole 
term in their gravity field. The non-uniform gravity field 
this allows is more physically appropriate for solid ob- 
jects (like TNOs) than the stars and giant planets for 
which the EKE01 model was developed. This non- 
uniformity is especially relevant for objects the size of 
most known TNBs (less than 400 km diameter), which 
are likely no t large enough to have reac hed hydrostatic 
equilibrium d Yasui and Arakawa , 2010l) . and thus could 
have relatively large quadrupole fields. The quadrupole 
moment is defined for an axisymmetric body as J2=(C- 
A)/(MR 2 ), where C is the moment of inertia about the 
polar radius, A is the moment of inertia about the equa- 
torial radius, R is the equatorial radius, and M is the 
mass. In the EKE01 model, the vector (X, Y, Z) pro- 
vides the angular precession rate relative to the iner- 
tial frame, and is in the (e m ,q mi h m ) orthonormal ba- 
sis, where e m is the normalization of the mutual orbit's 
Laplace-Runge-Lenz vector, which points in the direc- 
tion of the periapse, /z,„ is in the direction of the orbit's 
angular momentum vector, and q ln - hj„ x e,„. This 
vector due to solar torques is given by Equations 10 to 
12 in EKE01, and their relation to the secular evolution 
of the orbit is summari z ed by Equations A6 to A8 in 



Fabrvckv and Trema ine (2007P- iRagoz zine (2009J) for- 



mulated the additional precession due to the primary's 
quadrupole field as being: 



X 



2 ' [dbil (1 



2 \ a in 



lft"le 



-e 2 ) 2 



(1 -e 2 ) 

v in/ 



2\2 



a? 



(11) 



(12) 



4 \ a/„ 



d-4) 2 



2Q?, -Q? - Q? 



Where the terms a,-„, e-,„, and n m are the semimajor 
axis, eccentricity, and mean motion of the mutual or- 
bit, and R\ is the radius of the primary. In addition, 
Qii is the projection of primary's spin angular veloc- 
ity vector onto the axis i. Since the binaries dealt with 
in this work are of similar size, this process was re- 
peated to account for the secondary's quadrupole field, 
and the two quadrupole precession vectors added on to 
the {X, Y,Z) vectors defined in the EKE01 model. The 
total {X, Y, Z) vector was then combined with the dissi- 
pative terms (our Equation |7]l to feed the full evolution 
equations (Equations 1 to 6 in EKE01). 
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2.3. Integration Methods 



3. Monte Carlo Simulations 



Since the EKE01 model defines the evolution of the 
system by a set of four related inhomogeneous vector 
differential equations, we needed a numerical integrator 
that could solve them rapidly and precisely on a 
modern computer. For this we based our integrator on 
a Burlisch-Stoer method, which combines a modified- 
midpoint integrator with an polynomial interpolation 
method to increase precision and control error. Since 
the system was conservative when the dissipative terms 
where close to zero, we used a fixed timestep of 1.1 
mutual orbital periods when |Vi| + IV2I < 10" 18 sec -1 
(equating to an approximate circularization timescale of 
3.2xl0 10 years). On the other hand, if the system were 
dissipative, we used the adaptive tim e step m anagement 
algorithm described in Press et al.l (12007b . setting a 
minimum timestep of 10 days (864000 sec). This 
algorithm estimated the total error for each step as the 
root-mean-square of the normalized error estimates 
for each component (<?,„, e\ n , h m , h m , Qi, Qq), and 
kept it below a tolerance of 10~ 13 for each timestep. 
Throughout the simulation, the program keeps track 
of the total angular momentum (spin plus orbit) in the 
direction of the heliocentric orbit's angular momentum 
vector. As described above, this term determines the 
magnitude of solar perturbations on the mutual orbit 
and should be precisely preserved over the entire 
integration. The small number of cases (generally those 
with non-zero values of Ji) which did not preserve 
angular momentum were rerun at sufficiently smaller 
tolerances that momentum was again conserved. Con- 
versely, if this tight tolerance proved to be numerically 
unstable, the tolerance was slowly increased until it was 
stable but still conservative. 

We ended each simulation when it reached either 
4.5 Ga (i.e. the maximum physically possible evolution 
time) or reached an eccentricity smaller than 10~ 4 . 
This value of minimum eccentricity was chosen for 
an end state because preliminary simulations down 
to 10~ 10 showed no further orbital evolution beyond 
circularization. In addition, the simulation would end 
prematurely if the periapse fell below the Roche limit 
(«1.26(/?i + ^2), which we consider an impact), the 
apoapse grew beyond the Hill radius of the system, or 
one of the objects reached a spin period faster than the 
breakup rotation rate. The condition for the latter case 
was a rotation rate greater than 2n ^jGpl{3n), and in 
practice was never reached in our simulations. 



To test the responses of TNBs to the KCTF model, 
we conducted a series of Monte Carlo simulations in 
which we created a sample set of 1000 synthetic TNBs 
with randomized mutual orbital elements and system 
masses. The heliocentric orbit, physical properties, and 
range of rotation rates were all kept constant for each 
set. We then evolved each system for 4.5 billion years, 
or until the system either circularized, impacted itself, 
became unbound, or spun to breakup. 

As common initial parameters, we set the heliocentric 
orbits to a semimajor axis of 45 AU and eccentricity 
of 0.05, representative of the cold-classical belt, which 
contains the highest fraction of TNBs. We then varied 
the system GM range from 0.02 to 0.20 km 3 /s 2 . This 
corresponds to a radius range from 33 to 71 km for an 
equal-mass, p = 1 .0 g/cm 3 system, appropriate for the 
lower range of detectable TNB systems. The semimajor 
axis of the mutual orbit (in the frame of the primary) 
was varied uniformly from 0. 1 % to 10% of the system's 
Hill radius. This range is inclusive of nearly all known 
TNB orbits, as well as published formation models. 
The mutual eccentricity was varied uniformly from 
10~ 4 to 0.9999. The orbits were all orientated randomly 
on the sky by first generating a random direc tion for the 
hi„ vector with the algorithm of Knop lll970h . A second 
random vector perpendicular to the first was then 
generated and used as the e% vector. These two vectors 
plus the semimajor axis, eccentricity, and system mass 
then fully describe the randomized orbit. This random 
orientation cor responds to dynamical simulations of 
binary capture (Ko minami et all 1201 lb . Likewise, the 
initial spin poles of the two objects were pointed at two 
different vectors also randomized on the sky. 

We considered two different densities, 0.5 and 1.0 
g/cm 3 , representative of the range of reported densities 
for b inary TNOs in this size range ( Stansberrv et al.l 
20121) . As noted above, we also considered both 
solid and rubble-pile assumptions for ki. For each 
simulation set, we replicated the runs for each of the 
four combinations of density and ki. For most of the 
simulations, we allowed the initial rotational periods of 
the objects to vary between 4-48 hours. This range is 
consistent with observed lightcurves for solitar y TNO s 
falling and Bernsteini 120061 iThirouin et all 12010b . 
To check the sensitivity of this range, we also ran simu- 
lations with 2-7 day rotation rates. We ran most of our 
simulations with 72=0, but we also ran sets at 72=0.01 
and 72=0.1. The former being comparable to Saturn's 
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satellite Phoebe, and the later to Hyperion, both of 
which are in the same size range as our simulations. We 
ran most simulations with 2=100, a typically canonical 
value for large solid objects. Ho wever, since smaller 
objects may have much smaller Q dGoldreich and Sari , 
2009; IZhang and Hamilton!, l2008h . we also performed a 
set of simulations with 2=10- Finally, in addition to the 
equal-mass simulations, we ran a set with a mass ratio 
of 10: 1 . Since most observed binaries have a brightness 



difference of less than one magnitude ( Nol l et al 



2008a), these two mass ratios are representative of the 
observed population. 

The parameters varied per each set of 1000 simu- 
lations were thus Q, J2, rotation rate, mass ratio, 
density, and ki. Table Q] lists the 24 sets considered, 
for a total of 24,000 simulations and approximately 
1 million CPU-hours. Because the initial conditions 
of these systems were distributed across reasonable 
ranges of a, <?, i, and obliquity, they do not necessarily 
represent the primordial population of TNBs. Rather, 
they are a superposition of all the possible initial states. 
KCTF can then act as a filter to check whether a certain 
range of initial orbits (and thus formation methods) 
corresponds to the observed population of TNBs. 

Figure 3 here 

Table 1 here 

4. Results 

We found that a significant fraction of our simula- 
tions resulted in the synthetic TNB systems evolving 
to very tight circular orbits; Table [T] lists the relative 
fractions for each simulation. These tight, circular 
orbits were at less than 1% of rmii and eccentricities 
smaller than 10~ 5 , meaning that they were entirely 
dominated by mutual tidal interactions. In addition to 
these highly-evolved systems, a number of systems 
evolved to orbits that were tighter but still eccentric. 
Figure [2] shows the time evolution of one set of sim- 
ulations and Figure [3] separates initial and final states 
for the circularized and elliptical cases of the same 
set. The conservation of H' is quite evident in these 
plots, as is its dependence on cos I and -\fa. Figures 
lUE] and [6] show the final states of a range of different 
simulation sets, varying in both the speed of body tides 
and physical shape. 

Much of the tightly-bound population is beyond 
current TNB detection methods, including the WFC3 



camera on the Hubble Space Telescope and the laser 
guide star adaptive optics system at the Keck Observa- 
tory. However, some of the known TNO orbits do fall 
within this population. The closest published full orbit 
is (7 9360) Sila-Numan, f ormerly 1997 CS29, at 0.35% 



rum dGrundv et aU 120121) . The published eccentricity 



is 0.02, but a fully circular orbit is only excluded at 1 .8 
cr confidence . In addition, (120347) Salacia-Actaea at 



0.23% rum (IStansberrv et all 120121) and the centaur 



(65489) Ceto-Phorcys at 0.46% r HiU dGrundv et al 



2007b both have extremely circular orbits. These three 
objects likely represent the inner edge of the potentially 
very numerous tight circular population, and should 
become less unique as observations improve. 

Our simulations also naturally produced a deficit 
of systems with a final semimajor axis smaller than 
1.26 x (Ri + R2), our approximate Roche limit. It is 
perfectly possible that these systems could survive, ei- 
ther as a contact binary or brea king apart and reforming 
in a m ore stable configuration (IJacobson and Scheeres 
201 lb . However, any of these cases are beyond the ca- 
pabilities of our model and would require further work. 
This inner edge corresponds to a ///' of approximately 
0.4, where J is the total orbit + spin angular momen- 
tum, /' = ^GM^R e ff, M, is the system mass, and 
R e /f is the radius of a sphere of mass M, and density 
equal to the average density of the two objects. ICanup 
(2005) showed that binaries with J /J' < 0.8 can be 
produced by collisions. Thus, we found that KCTF can 
transform a wide range of binaries formed by capture 
into ones sufficiently close as to be indistinguishable 
from collision-produced binaries, t he general case o f 
what was shown for Orcus-Vanth by RagozzineJ (2009). 

Because the quadrupole component of Kozai per- 
turbations are axisymmetric, the resulting systems also 
preserve the direction of their initial mutual orbit rela- 
tive to the plane of their heliocentric orbit. The EKE01 
model only includes the quadrupole component, though 
in the actual dynamics of these systems there are also 
higher-order terms. The next higher term is the octupole 
component, which is not axisymmetric and can there- 
fore cause a syste m to flip from pr ograde to retrograde 



(and vice versa, Naoz et al.J, 1201 II) . However, the 



relative strength of the octupole to the quadrupole 
goes as ai n /aheiio, and so the quadrupole completely 
domi nates in all the cases we considered dRagozzinel 
2009). Thus, the initial prograde/retrograde ratio was 
preserved in all our simulations. 



6 



Figure 4 here 

4.1. Stability of Orbits to KCTF 

A further inclination effect can be seen in Figure |3j 
the inclination region within 10° of perpendicular to 
the heliocentric orbit is empty for all but the tightest 
semimajor axes. This range of inclination equates to 
such low values of H' that all orbits starting there reach 
very high levels of eccentricity, and therefore initiate 
runaway KCTF decay (e.g. Figure[TJ. Figure [5] clearly 
shows that these orbits mostly end up tighter than 1 % of 
rmu- Indeed, as Figure|5]shows, there is a clear limiting 
value of \H'\ below which orbits are not stable to KCTF. 

We define this limiting value, \H' tlde \, as the mini- 
mum value of \H'\ where the relative difference in the 
initial and final semimajor axes is less than 10%, i.e. 
\cij - a/\/ai < 0.1. The values of \H' dd J for all our simu- 
lations are listed in Table Q] The average values of this 
\H' lide \ are all close to 1.0, with the higher-dissipation 
cases above and the lower below. Thus, \H'\ appears to 
be a good normalized indicator of KCTF susceptibility 
for TNB orbits; values below 1.0 likely will have 
experienced Kozai cycle-driven tidal evolution, while 
values above will likely have not. This is especially 
useful for determining the stability of observed binary 
systems where the mutual orbit is known, but not any 
physical or rotational properties. 

In addition, we found it useful to define \H'. \, 
the minimum value of \H'\ for an orbit that did not 
circularize. This value ranges from 0.2-0.7 for non-^ 
runs, but is close to zero for the simulations with /2- 
This shows that Ji can provide an island of stability 
for eccentric orbits close enough that Ji blocks fur- 
ther Kozai cycles, but not close enough for further 
tidal decay. These systems have tidally evolved before- 
hand to reach this state, but are stable once they reach it. 

Figure 5 here 

4.2. Effects of Physical Parameters 

For our base sets of simulations, we assumed 2=100, 
72=0, rotation rates between 4 and 48 hours, and equal 
masses for the two objects. Using these parameters, we 
ran sets at each combination of density equals 0.5 or 
1.0 g/cm 3 and a Love number appropriate for either a 
rubble pile or elastic solid body. The remainder of the 
simulations perturbed one of the first set of parameters 
and ran the same set of four density/Love number 
combinations; Table[T]shows this grouping. 



In general, the low-density rubble piles were the 
most susceptible to tidal decay, followed by the high- 
density rubble piles. This makes sense, as the internal 
friction of a rubble pile allows it to dissipate tidal 
forces very effectively. The lower density allowed a 
larger radius per mass, thus raising larger tides and a 
higher Love number (by Equation @J. Less obvious is 
that fewer of the low-density elastic solid simulations 
decayed to stable circular orbits than the high-density 
cases. Here, the difference can be seen in the final 
column of Table [1] these low-density, rigid systems 
suffered a much higher rate of mutual collisions. Their 
larger radius prevented a significant fraction of the 
very eccentric systems from circularizing, lowering the 
overall efficiently of producing close, circular orbits. 
As noted below, these mutual collisions could still 
produce close or contact binaries, but that is beyond the 
capabilities of these simulations. 

While we assumed the canonical value of tidal 
2=100 for most of the simulations, a lower value is 
likely more physical for t he si ze of objects we consid- 
ered (IGoldreich and Sari, 2009). The simulations that 
we ran with 2=10 did show an average of 20% greater 
propensity to tidally decay and circularize. This small 
increase in tidal decay for a full order of magnitude 
decrease in Q shows the true driver in this evolution 
is the closest periapse the system reaches. Either the 
objects become close enough to undergo decay or they 
do not, there is not much space in between. Indeed, the 
increase is roughly comparable (lOO/lO) 1 ^ 8 , where the 
strength of the tides goes as a 8 . An additional effect of 
the lower Q is to remove the higher mutual collision 
probability for the low-density, elastic solid case. As 
shown in Figure |4] the faster tides allow these cases 
to circularize at larger separations, limiting the impact 
probability. 

We ran sets of simulations with a Jo of both 0.01 
and 0.1. As seen in Figure [6] the main effect of Ji 
on wider systems was to decrease the obliquity of the 
objects with respect to their mutual orbit. For slightly 
more evolved systems, though, an interesting effect 
could be seen; a number of the systems decayed to 
a point where J2 was strong enough to block further 
Kozai cycles, but their periapse was wide enough at 
the point of being frozen that no further tidal decay 
could occur. This freezing-in of Kozai-blocked systems 
created the island of high-inclination, moderately 
eccentric systems seen in Figures [3] and [5] Since the 
frozen systems did not evolve further, they reduced the 
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total fraction of systems that circularized. However, 
since they had to first tidally evolve to reach the frozen 
state, the total fraction that tidally evolved (and thus 
\H'- |) is almost the same as the non-72- Also, since 
tidal evolution is strongly driven by the periapse, there 
is a noticeable gradient in eccentricity, with the closest 
frozen systems having e < 0.2. Thus, the effect of J2 
on KCTF of TNBs is produce stable high-inclination, 
moderately eccentric TNBs, while not eliminating the 
tight circular population. 

We also ran a group of simulations with much 
slower initial rotation rates, from 2-7 days (and at 
random initial orientations, just as before). These 
simulations were slightly less likely to evolve to close 
circular orbits than their fast rotating counterparts. This 
may be because the final orbital periods of the close 
circular orbits were closer to the initial rotation rates of 
the faster rotators. Thus, the slower rotators required, 
on average, a larger amount of spin-orbit interaction to 
reach a doubly-synchronous state (with both spin peri- 
ods equal to the orbital period). The effect was slight 
enough, though, to conclude that initial rotation rate 
is not a significant factor in the evolution of these TNBs. 

Finally, we changed the mass ratio of the two ob- 
jects to be 10:1. Since most known TNBs have 
near-equal brightness ratios (tNoll et al.i |2008a), this 



covers most of known TNBs that are not dwarf plan- 
etary systems. The effect on tidal decay was again 
similar to the slower initial rotation rates, a slight 
decrease in the efficiency of producing close circular 
orbits. Here the reason is that most of the dissipation 
is taking place on the secondary, which can only have 
smaller tides due to its smaller radius. Again, though, 
the effect is slight enough that most systems would be 
completely insensitive to mass ratios between 1:1 and 
10:1. 

Table 2 here 

4.3. Survival Probability 

Not all of the simulated binaries survived to either 
circularize or reach 4.5 Ga. As shown in Table [1] about 
1-15% of our random initial orbits proved in some way 
unstable. The largest fraction of destroyed systems was 
due to impacts when Kozai Cycles drove the periapse of 
the system within the mutual Roche radii of the objects. 
This was especially apparent for systems with elastic 
solid values of Icl (weak tides) and large radii, which 
consistently had destruction rates higher than 10%. Ev- 
ery other case had destruction rates smaller than 10%. 



As noted above, this is because the low-density solid 
systems circularize slowly enough that semimajor axis 
decay can bring their periapse below the Roche limit. 
In real TNBs, such a close approach could result in the 
temporary breakup of the secondary and reformation of 
one or more new satellites in tight circular orbits. The 
simulated systems were also considered destroyed if 
the apoapse of the system exceeded ruuu though this 
case did not occur in any of our simulations. Likewise, 
we had no systems spin to breakup, with spins either 
slowing to synchronize or not changing at all (see 
Figure |6). It therefore appears that mutual collisions 
are the only practical way to destroy a TNB system 
with KCTF, and even then only with an efficiency of a 
few percent. In addition, 90% is a lower limit for the 
number of TNBs (of an initially random population) 
which survive as binaries once formed, exclusive of 
impacts or close encounters. 

The survival rate i n our simulations i s much higher than 
that calculated by iPetit and Mousisl (120041) . who found 
that non-disruptive impacts were the most effective 
means of destabilizing wide TNBs. This implies that 
small impacts and disruptio n by Neptune encounter 
( Parker and Kavelaars . 2010h are still the most efficient 
means to cause a TNB to become unbound. However, 
our results show that KCTF is very efficient at trans- 
forming wide binaries into tight ones, which would 
presently be observed as single objects. Indeed, KCTF 
can quickly shrink wide binaries down to orbits tight 
enough to have a high su rvival probability in the even t 
of a Neptune encounter (Pa rker and Kavelaars , 2010l) . 
Thus, any scaling of the initial population of wide 
TNBs based on survival rates must also account for the 
fraction which decay into much tighter orbits. 

Figure 6 here 

4.4. Obliquity 

All the systems we simulated started out with 
randomly directed spin axes, and thus each object had 
random initial obliquities. We define obliquity here 
to be the angle between the spin pole of the objects 
and a vector perpendicular to the plane of the mutual 
orbit. Thus, an obliquity of 0° is parallel to the orbit's 
angular momentum vector, while an obliquity of 90° is 
perpendicular. Figure [6] shows the final obliquities of 
equal-mass, 2=100 simulations with and without Ji. 

In general, the non-72 runs experienced two stages of 
obliquity evolution. At less than about 3% rnm, the 
spins of the two objects begin interacting with each 
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other and with the orbit. This causes any objects rotat- 
ing retrograde to their orbit to flip around to prograde, 
causing the step jump in Figure [6] In addition, the two 
objects begin to match their obliquities and rotation 
rates to each other, thus producing the star shapes in 
Figure [6] The second (and terminal) stage is for the 
objects to match their spin poles and rotation rates to 
their mutual orbit. This drive to zero obliquity only 
occurs for the non-/2 simulations when tidal forces are 
also strong enough to circularize the orbit. 

The simulations with J2 were similar, but the ad- 
ditional gravitational factor from J2 caused them to go 
to zero obliquity much faster. Indeed, even some very 
wide orbits with little tidal evolution were still driven 
to zero obliquity by J2- Shape effects thus can have a 
major effect on the rotational evolution of even wide 
systems. 



A further effect can be seen in Figure [6] as some 
circularized orbits ended with non-zero obliquities and 
non-synchronous rotation rates. This effect was more 
pronounced for less-dissipative simulations. These 
systems appear to have been captured into Cassini state 
2, which is low obliquity for f ast precession and hig h 
obliquity for slow precession dFabrvckv et all 2007 ). 
As they circularize at near-constant semimajor axes, 
their precession rate due to Jo drops (equations 1 1-13). 
Since they are sufficiently close that J2 precession 
dominates, the overall precession rate becomes very 
small and the objects follow Cassini state 2 to high 
obliquity. 



This capture process is not necessarily observable 
in real systems. Most small objects with known shapes 
(asteroids and minor satellites) are not pure oblate 
spheroids, but rather more complex shapes. It is 
possib l e for these shapes to capture int o Cassini state 2 
ilPealel 1 19691: iBills and Nimmol 120081) . but it is much 
more difficult. In addition, they may have precession 
rates due to higher-order terms that are large enough 
to prevent even a very circular system in Cassini state 
2 from going to high obliquity. Thus, a much more 
through study of the possible post-circularization 
rotational evolution of TNBs is needed to properly 
predict their current states. 

Table 3 here 



5. Discussion 

KCTF provides an evolutionary path to convert wide, 
elliptical binaries into close, circular ones. It therefore 
helps to explain the observed dichotomy between 
the few (but well-sampled) wide TNB systems and 
the apparently numerous tight systems. In addition, 
because it is a process that does not require any external 
forces other than the Sun, it applies to all objects that 
could be called trans-neptunian binaries, from classical 
Kuiper Belt objects to highly scattered centaurs like 
Ceto-Phorcys. And, because it is independent of the 
surrounding disk, KCTF can shrink and circularize or- 
bits over a much longer timescale than disk dynamical 
friction (e.g. L 2 S). 

A consequence of this extensive evolution is that 
KCTF should have significantly reshaped the orbits of 
most TNBs since their formation. In the process, some 
information is lost, as tidal decay is an irreversible 
thermodynamic reaction. The distribution of semimajor 
axes for present-day TNBs is thus necessarily tighter 
than at formation. By what factor it is tighter is hard to 
determine, as the physical properties of the objects do 
affect the efficiency of semimajor axis decay (see Table 
[TJ. Similarly, the exact eccentricity and inclination of 
the orbit are changed in ways that have partially erased 
their history. 

More clear is that KCTF generally preserves the 
prograde/retrograde ratio of initial distribution. Though 
the octupole term from th e solar per t urbati ons can 
flip a very inclined system dNaoz et all 1201 lb . it also 
has no preference fo r the di rectionality of the system. 
ISchlichting and Sari (2008) predicted a retrograde 
preference for binaries formed by the L 2 S capture 
method. On the other han d , the g ravitational collapse 
method of INesvornv et al.l (120 101) favors orbits in the 



direction of disk clump rotation, which simulations 
(cited therein) generally show as prograde. Three-body 
methods (L 3 and momentum exchange) have a much 
weaker inclination dependence. The current inclination 
distribution could therefore serve as a tracer for these 
various formation methods. Because the sense of 
motion can be hard to determine for TNBs, only a 
few TNBs have published unambiguous inclinations. 
However, for the 16 TNB orbits listed in Table [2] there 
is a 4: 1 prograde preference. 



The H' parameter introduced above is also pre- 
served by KCTF (except for spin-orbit interactions), 
and Table|2]also lists the H' for those 16 known orbits. 
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The widest orbits (a > 5%rHiii) are apparently stable, 
with their large a and low inclinations keeping them 
stable even at high eccentricities. These systems likely 
have thus been undergoing low-amplitude Kozai cycles 
for most of the history of the solar system. Table |3]lists 
the known TNB systems in near-circular orbits; these 
systems appear to represent the outer edge of the close 
circular population identified in our simulations. As 
more sensitive high-resolution imaging systems come 
on line, an increasing number of these tight systems 
should be detected. The remaining systems smaller 
than 5%rnm have values of \H'\ at or below unity, 
implying that they all should have some amount of tidal 
evolution. Most still have a larger \H'\ than the typical 
values of \H'- \ in Table [T] meaning they could have a 
range of physical properties and still not circularize. 



The inclusion of Ji lowers the effectiveness of 
KCTF, but does not eliminate it, especially for rubble- 
pile objects. In addition, creates an island of stability 
that allows otherwise unstable observed system to 
be in permanent eccentric orbits. A slower initial 
rotation rate or 10:1 mass ratio also slightly lower the 
effectiveness of KCTF, but do not change the basic 
trends. 

The observed population of TNB orbits fits well 
to our simulations with J 2- These simulations predict 
that, as high-resolution observational systems improve, 
a large number of TNBs will be detected with very 
tight, circular orbits. Indeed, considering the fraction 
of known wider binaries, tight near equal-mass TNBs 
may be extremely common. 



Two eccentric systems do have an \H'\ small enough 
to stand out, though: 2004 PBiog and 2001 QC 29 8- 
These two systems happen to be the only scattered-disk 



binaries by the DES classification dElliot et al.L 2005 



Osip et all 120021) in Table [2] but we think it more likely 
that their apparent stability is due to their non spherical 
shapes. Indeed, TableQ]shows that if both of the objects 
in each system had a J2 of at least 0.01, they could 
easily be stably eccentric for the lifetime of the solar 
system. The other eccentric systems do not require 
J2 to be stable, but it would do no harm to their stability. 

The three systems just wider, 2001 XR254, Altjira, 
and 275809, all have values of \H'\ just larger than one. 
Since they are all also in the classical Kuiper Belt, these 
systems have potentially had very little tidal evolution 
since their original formation as binaries. Interestingly, 
they also have very similar systems masses and mutual 
orbit inclinations. Future investigations could help 
determine if this is simply coincidence or an optimal 
point for stability. 

6. Conclusions 

KCTF can significantly transform the orbits of trans- 
neptunian binaries. At least 90% of random synthetic 
TNB systems survive 4.5 Ga of KCTF evolution. A 
third to half of the surviving TNB systems decay to 
circular orbits at less than 1% of their mutual Hill 
radius. Some of these systems can have values of J/J' 
similar to impact-generated systems. The remaining 
systems are stable being eccentric over the lifetime of 
the solar system. All resulting systems preserve their 
initial prograde/retrograde preference. 
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Table 1: Summary of KCTF simulation results: % Circ. are the fraction of surviving orbits with e < 10 , % Dest. are the fraction of initial orbits 
destroyed over the simulation, and IffijJ and \H'.\ are as defined in the text. 
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Table 2: H' for observed sys tems with fully-cons t rained orbits. / is the angle between the heliocentric and mutual orbital planes, and GM is the 
system mass. Orbits are from Grundv et al. (2011), Parker et al. (2011), Shepnard et al. (2012) and references therein. 
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H' 
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0.35 


0.02 
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0.72 


-0.32 


2001 QC 2 98 




0.50 


0.34 


73.5 


0.79 


+0.19 


66652 


Borasisi 
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49.4 


0.23 


+0.55 
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+0.84 


58534 


Logos 


3.25 
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0.16 


-1.43 
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Table 3: Observed TN B systems with near-circular orbits. GM is the system mass. Orbits are from iBrown et al l ( 120101) . iGrundv et all 1201 lh . 
Stansberrv et al. (2012) and references therein. 



Designation 


Name 


a (% r m u) 
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GM (km 3 /s 2 ) 


120347 


Salacia 


0.23 


0.01 
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0.02 


0.72 
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Figure 1 : An example of rapid KCTF evolution of a TNB mutual orbit; see Section 2. 1 for a full description. 
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Figure 2: The evolution over time of 1000 synthetic equal-mass TNBs with Q = 100, p = 1.0 g/cm , elastic solid and 72 = 0. 
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Figure 3: A comparison of end states for three different physical properties and observed orbits. The simulations are all equal-mass, Q 
Observed orbits are as listed in Table [2] 
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Figure 4: The circularization times for three different physical properties; the simulations are all equal-mass systems. The horizontal line corre- 
sponds to 4.5 billion years, and any points along it represent stable elliptical orbits. The points on the left below the line are simulations that were 
evolved by KCTF to close circular orbits, while the points on the right below the line started at very low eccentricity and evolved to e < 10~ 4 . The 
final semimajor axis in tidally-evolved systems is usually very close to the periapse of the orbit when body tides become the dominant force, as 
most of the energy loss in an evolving orbit is at the periapse. 
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Figure 5: A comparison of initial \H'\ and final semimajor axis for three different physical properties and observed orbits. The simulations are all 
equal-mass, Q = 100. Observed orbits are as listed in Table [2] 
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Figure 6: The final spin rate and obliquity (to the mutual orbit) with and without Ji ; upward-triangles are the primary objects, and downward are 
the secondaries. The high obliquity of some of the close orbits with Ji is due to them being captured into Cassini state 2 (see text). 
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